#!/bin/bash

function delFile {
	if test -s $1; then \rm $1; fi
}

function getIsoformWithNewJuctions {
	if [ $# -ne 4 ]; then
		echo "This shell get the junction fragment seq and check if a *.fa"
		echo "file containing the seq of fragment. It needs arg"
		echo "(1) the sequence (.fa)"
		echo "(2) the junctions file"
		echo "(3) the isoform file"	
		echo "(4) the output file"
		echo "4 ptrs are needed"
		echo "$# is provided" 
		echo "$@" 
		return 0
	fi
	local faFileN=$1
	local juncFileN=$2
	local shortFragFileN=$2'.seq'
	local isoformFileN=$3
	local newIsoformFileN=$4
	if test ! -s $faFileN; then
		echo "$faFileN doesn't exist"
	elif test ! -s $juncFileN; then		
		echo "$juncFileN doesn't exist"
	else
		python ~/SCRIPT/get_junction_seqs.py $faFileN $juncFileN $shortFragFileN 
		#>> /dev/null
		if test -s $isoformFileN; then
			python ~/SCRIPT/get_seqs_w_short_fragments.py\
			 $isoformFileN $shortFragFileN $newIsoformFileN 
		#>> /dev/null
			\rm $shortFragFileN		
		fi
	fi
}

# for test 
#gene=LOC728661
#new_junc_file='1_'$gene'.bed'_'2_'$gene'.junctions_.except'
#getIsoformWithNewJuctions $PWD/tmp/$gene'.fa' $new_junc_file $gene'_I_selected.fa' 'newPep.fa'

# For each gene list, concatenate the table with each column as
# chr gene foundJunctionNo candidateIsoformNo selectedIsoformNo (total 5 columns)

function  get_new_peptide
{   	    				
	while read line
	do
	    if [[ $(expr match $geneList "#") == 1 ]]; then
	        continue   # comment 
	    fi
	    # echo will turn tab into space automatically
	    # echo $line
	    local gene=$(echo $line | cut -f 1 -d' ')
	    local chr=$(echo $line | cut -f 2 -d' ')
	    if [ -d $gene ]; then    	
	    	cd $gene
	    else 
	    	echo "$PWD $gene doesn't exist "    	
	    	continue
	   	fi
	    # TODO have a better naming of the junction file
	    local new_junc_file='1_'$gene'.bed'_'2_'$gene'.junctions_.except'	      	 	
	    if ! test -s $new_junc_file; then
	        cd ..
	        continue # no junction are found 
	    fi
	   	local faFileN=$PWD/tmp/$gene'.fa'
		local juncFileN=$new_junc_file		
		local isoformFileN=$gene'_I_selected.fa'
		local newIsoformFileN=$chr_$gene'.tmp'
	    getIsoformWithNewJuctions $faFileN $juncFileN $isoformFileN $newIsoformFileN
	    if test -s $newIsoformFileN; then	    	
	    	echo -n ">"$chr"_"$gene"_" >> $outFile
	    	cat $newIsoformFileN >> $outFile
	    	\rm $newIsoformFileN	    	
	    fi
	    cd ..    
	done < $geneList	
}

if  test ! -s "$1" ; then		
    echo "This shell will extract the isoforms which consist of the new junctions"     
    echo "The 1st prt is the file name of a gene list."
    echo "The 2nd prt is the output file"    
else 
	geneList=$1
	outFile=$2
	get_new_peptide $geneList $outFile
	noNewPep=$(grep '>' $outFile | wc -l | cut -f 1 -d ' ')
	echo "$geneList has $noNewPep new peps."
fi 
	
	
	
